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Abstract. We present a model for the lower overshoot layer of the Sun, based on the realistic solar stratification, 
without the use of a "m ixing-length" paramet e r, by solvin g the system of Reynolds momentum equations using 
the closure formalism of ICanuto fc Dubovikovl il997l Fl998'l. A fixed value of velocity anisotropy is assumed, and 
the local convection model is assumed to be valid for the convectively unstable layer. In accordance with seismic 
constraints, overshoot (defined as the amount by which the convectively mixed zone extends beyond its boundary 
in local theory) is found to be as low as about 6 percent of the pressure scale height, and it is not bounded by a 
discontinuity from below. 



1. Introduction 

The layers situated below the convective zone of the Sun 
are currently considered to play an important role in the 
dynamo m echanism that m aintains the general solar mag- 
netic field l|PetrovavL l2000). Besides, turbulent mixing pro- 
cesses taking place in this layer are invoked to explain 
the lithium deficit of the solar atmosphere and to explain 
the remaining small disc repancies between t he seismic and 
standard solar models l|Gough et all fl996'). This turbu- 
lence may result partly from convective overshoot, and 
partly from rotational instabilities in the solar tachocline. 
The relative importance of these two sources of turbulence 
is an important question. 

The first models for the lower overshoot layer were 
based on nonlocal extensions of mixing-length t heory 
llvan Ballegooiienlll982tlSkalev fc StixUl99lHZahrlll99lt 
lUmezuL 11992'!. These predicted 1 that the thermal strati- 
fication of the overshoot layer is nearly adiabatic, its ex- 
tent is a significant fraction of the pressure scale height 
Hp, and it is bounded by a discontinuity in the superadi- 
abaticity AV = V - V ad from below. (V = d\nT/d\nP, 
and V a( j = (dlnT/dln P)s, as usual, P, T and S being 
pressure, temperature and entropy.) 

Such a discontinuity should, however, lead to a no- 
ticeable increase in the amplitude of the oscillatory sig- 
nal i n the frequency distribution of the solar oscilla- 
tions ijChristensen-Dalsgaardl Il99fj|). Attempts to detect 
this have remained unsuccessful ijMonteiro et all Il994t 
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In fact, adiabatic stratification is explicitly assumed in 
some nonlocal MLT models, while in others it is a direct con- 
sequence of the assumption of a perfect correlation between 
vertical velocity and temperature fluctuations. 



IChristensen-Dalsgaard et all Il99st iBasu fc AiitiaL l200l|) . 

which implies that the non-local mixing-length models are 
incorrect. Thus, either the overshoot layer is not bounded 
by a discontinuity from below , or its full extent is very 
small, according to lBasul l)l997l) not more than 5 % of Hp. 
It should be stressed, however, that even in the case of 
a discontinuity, this limit is only valid for overshoot in 
the helioseismic sense, i.e. for the amount by which con- 
vection extends beyond its boundary calculated in a local 
model. 

Note also that, as ICanutol l|l999t) pointed out, from 
the point of view of stellar evolution the kinematical or 
dynamical extent of the overshoot, i.e. the extent of the 
mixed layer, is a more important quantity than its ther- 
mal extent, even if the latter were to be adiabatic. Thus 
throughout this paper we will define the lower boundary 
of the overshoot layer at the point where the kinetic en- 
ergy is reduced by two orders magnitude compared to its 
value at AV = 0. 

A number of numerical simulations have also been car- 
ried out to study the problem, but they have not been con- 
clusive. Owing to the extremely different dynamical and 
thermal timescales in the solar interior, the simulations 
are eit her very far from the real solar interior in parameter 
space ( Roxburgh & Simmon slll99fll Hurl burt et al 
ISingh et allll995HSaikia et alll2000l) or they never reach 
a the rmally relaxed state l)Andersenl ll994l iFrevtag et all 
Il99t1) . 

In these circumstances, the most promising way to 
model the overshooting layer is the Reynolds stress 
method, based on the well known Reynolds momentum 
hierarchy of the hydrodynamical equations for a turbu- 
lent medium. Our aim in such a model is to determine 
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the radial distribution of the turbulent kinetic energy fc, 
the mean square relative temperature fluctuation q, the 
normalized energy flux J and the energy dissipation rate 
e under the convectively unstable region. Note that a de- 
termination of e is equivalent to the determination of the 
mixing length I as this parameter is usually introduced 
by defining / = fc 3 / 2 je. Thus, in this approach the mixing 
length is determined consistently within the model and it 
is not treated as a free parameter. 

The Reynolds moment equations for fc, q and J were 
solved for the lower oversho ot la yer by lUnno fc Kondol 
dl989h . lXiong fc Chenl lll992l) . and lXiong fc Penal fcOOllT 
resulting in significant non-adiabatic overshoot. However, 
in those calculations / was treated as a free parameter, 
arbitrarily set to be equal to Hp. The use of an equation 
for e to get r id of the free parameter was suggested by 
ICanutol (l993) and it has b een applied in a simplified k— 
e model bv IPetrovavl l|l998j) . A fully consistent Reynolds 
stress formalism for stellar convection has been developed 
in a series of papers by Canuto and co-workers. In this 
paper we will employ a formalism basica lly identical to 
that of ICanuto fc Dubovikovl l|l997l ll99S|) . This formal- 
ism has been shown to agree well with numerical simula- 
tion results for a generic convective flow with ineffi cient 
convection and unit Prandtl number l|Kupkal . |l999^ . 

By solving the momentum equations for the realistic 
solar stratification, here we present the first model of the 
lower overshoot layer without an assumed mixing length 
parameter. 

2. Equations and closure 

We will use the following notations. Any variable / is split 
into a mean and a fluctuating part as / = / + /'. The 
velocity v has only a fluctuating part. We assume plane 
parallel geometry and the depth z is measured from AV = 
0, so the gravity acceleration g is positive. Beside g for 
density, we also introduce the notations q = (T'/T) 2 , J = Q z p 



ICanutol fl993t iGrossman et all Il993|) . With the assump- 
tions and notations mentioned above, and after customary 
dimensional modelling of some terms, they read 



d t k = -d z F k + S p gJ - e 



dt-J = -d z Fj - C qJ S p gq 
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where the non-local fluxes are 



F k = wv 2 /2 



F a 



, - W {T'/Tf 

Fj = w 2 T'/T, F e =wTi 



(2) 
(3) 

(4) 
(5) 



£i being the local dissipation rate. For the coefficients we 
use the following values: C k .j = 1.0, C q j — 1, f a = 1, 
C e = 1.92, C j£ = 1.44, 5 P = 1. 

These e quations are essentially a special c ase of the 
equations of ICanuto fc Dubovikovl l|l997l Il998|) for a thin 
layer, and, for simplicity, with fixed anisotropy. (Note, 
however, that our definiton of J and q slightly differs from 
theirs.) Correspondingly we supplement Eqs . Q -ffl with 
the closure suggested bv lCanuto fc Dubovikovl l|l997|) : 
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w(T'/T), w — v z , k — v 2 /2 — w 2 j f a where f a is an 
anisotropy parameter. 

We make the following assumptions: 

— e7e< 1 

— The Reynolds number Re 3> 1, while the Prandtl num- 
ber y/x< 1 

— The extent of the overshoot layer d <C Hp , where Hp 
is the pressure scale height. 

— The turbulent flow field is characterized by a mild and 
fixed value of the anisotropy f a S 0(1). This assump- 
tion is made for simplicity only. 

With these assumptions we can use the Boussinesq ap- 
proximation q 1 I q — —5p(T'/T) where Sp is an order of 
unity factor (<5p = 1 for full ionization). Following stan- 
dard practice, we will also neglect the vVP' term in the 
equation for fc. 

The detailed derivation of the Reynol ds momentum 
equations was given by many authors IjXiond . Il980t F c + F k + F r 
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where A = K 2 e- 4/3 , v t = ^K 2 /e, a t = 0.72 and \ is the 

radiative conductivity. Since wT' = JT and T' 2 = qT 2 , 
we can obtain d z Fj and d z F q multiplying d z Fj> by 1/T 
and d z F 1/2T jj by 2/T 2 . 

Finally, the t imescales t^h and t„ were calc ulated using 
formulae (34) of lCanuto fc Dubovikovl l)l998|) . 

3. Numerical solution 

We are looking for stationary solutions of Eqs. CQ)-lQ}, in 
which case the system reduces to an eighth-order system 
of ordinary differential equations. The condition of flux 
equilibrium 
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where Fq is the total solar flux, F c is the convective heat 
flux, Fk is the kinetic energy flux, and F r is the radiative 
flux, must also be coupled to the problem. 

The upper boundary conditions are set in a point 
slightly within the convectively unstable layer, where we 
fit the solution to the solution given by a local convection 
theory. The equations of local theory are identical to Eqs. 

except that the non-local fluxes are not present 
and the radiative terms are negligible. So these equations 
reduce to an algebraic system for our variables. As the 
system is redundant, in addition we also use the condition 
HI 0|l to set the boundary values. The only free parameter 
is AV which was given a small positive value (~ 10 -8 ). 

A relaxation method was used to solve this nonlinear, 
coupled system of differential equations. We start the in- 
tegration from an arbitrary initial distribution where the 
dependent variables are gaussian functions of z, of a width 
comparable to Hp. During integration the following re- 
strictions were imposed: 

— Below the point where fc = or e — 0, all dependent 
variables are set to zero. 

— Wherever division by fc or e occurs in the equations we 
have used the condition k > 10" 100 and e > 10" 100 , 
to avoid floating point overflow. 

— At points where the velocity-temperature correlation 
J = J/ifakq) 1 / 2 > 1 or < —1, this quantity is set 
to 1 or —1, respectively. Note that this only affects 
the (physically irrelevant) intermediate stages during 
relaxation, while in the final converged solution, the 
unphysical situation mentioned above does not occur. 

In the initial state we use a linear profile for AV 
(AV = — 0.1 z/Hp); this crudely represents the actual run 
of the thermal stratification below the convectively unsta- 
ble layer. The iteration goes by making one "timestep" 
using the full, time-dependent set of Eqs. 
lowed by an instant adaptation of the temperature pro- 
file to ensure the condition of the constancy of the solar 
flux. This latter is achieved by calculating the kinetic flux 

(Ff. — 1/2 gw 2 ) and the convective flux (F c — C p gwT'), 
and subtracting them from the known solar flux to get the 
radiative flux 

7 3qk dr ' 

from which we can calculate a new temperature stratifi- 
cation (k is the opacity and a is the Stefan-Boltzmann 
constant). 

The whole process is repeated until it converges to 
a stationary state. The relaxation does not describe a 
real temporal evolution, as the instant adaptation of the 
temperature profile is unphysical, but in this manner we 
can circumvent the problem of forbiddingly large radia- 
tive/dynamical timescale ratios. Thus, the technique is 
essentially an artificially accelerated thermal relaxation. 

As the layer is thin, the perturbation of the hydrostatic 
equilibrium is expected to remain small. The background 
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Fig. 1. Distribution of turbulent kinetic energy fc, relative 
squared temperature fluctuation q, convective flux equiva- 
lent speed J, and energy dissipation rate e in our model, as 
functions of depth below the convectively unstable layer. 
(All dimensional variables here and in other figures are 
given in CGS units.) 

values of g(z), p(z), k(z) and Hp(z) were thus assumed to 
be unaffected by overshooting, and they were taken from 
a mor e recent version of the solar model of ldienther et all 
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Fig. 2. Distribution of the turbulent length scale ("mix- 
ing length") I, time scale r, and superadiabaticity AV in 
our model, as functions of depth below the convectively 
unstable layer 



4. Results and discussion 

Our stationary solution is presented in Figs. 1 and 2. The 
turbulent length and time scales in Fig. 2 are defined as 
I = fc 3//2 /e and r = fc/e. The thickness of the overshoot 
layer is strikingly small, below six percent of Hp, and the 
AV curve does not show a particularly sharp disconti- 
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nuity. Th is is compatible with the available helioseismic 
evidence llBasu fe AntialEooH) . It must be stressed that, 
by design, our models can only yield the extent of over- 
shoot in the helioseismic sense (i.e. the amount by which 
the convective zone extends beyond its boundary in local 
theory). It is quite possible that a proper nonlocal treat- 
ment of the full convective zone would lead to slightly 
subadiabatic stratification and to a downward directed F c 
(the more conventional criteria for overshoot) already well 
above the boundary of the convective zone as defined in 
local theory. In this sense, the total extent of the over- 
shoot layer may far exceed our estimate of the seismic 
overshoot. 

It is interesting to note that the "mixing length" I 
is found to decrease continuously towards the bottom of 
the overshoot layer. This underlines the incorrectness of 
the assumption I — const, in those models that use a free 
length parameter. 

A doubt may arise regarding the validity of these find- 
ings, given that our model only extends to the overshoot 
layer itself, assuming the validity of the local description 
for the unstable layers above. One sign of the artificiality of 
this sudden introduction of nonlocality is the sharp jump 
in the J curve from negative to positive values right after 
z = (not visible in the figure). This "patching" of our 
nonlocal overshoot model to the model of the SCZ com- 
puted with a local convection theory also implies that e.g. 
the full extent of the convective zone cannot be determined 
independently from our model. Indeed, placing the top 
boundary of the computational domain to a higher level 
(z = —Q.2Hp) we found that overshoot in this case started 
right at the upper boundary again. The total extent of this 
"seismic" overshoot, however, showed only moderate sen- 
sitivity to the placement of the upper boundary (resulting 
in an overshoot distance of 0.06-ffp in the case quoted 
above). A model incorporating the full convective zone, 
with a nonlocal treatment throughout, is thus clearly a 
logical future extension of the present model. However, the 
solution presented Fig. 3 may give us a reassurance that 
the small overshoot is not an artefact of our method. In 
this plot the stationary solution of Eqs. JU-JU is shown, 
obtained using the simple down-gradient closure 

F k = - (k 2 /s d z K) Fj = - (k 2 /e d z j) (11) 
F g = - (k 2 /e 8 z q) F s = - (k 2 /e d z e) (12) 

instead of Eqs. ©-JHl)- Note that this closure is known to 
be wrong, and it is only used here for comparison. It is 
evident that a deep overshoot can also be achieved by our 
method of solution; thus, the nature of the closure is the 
chief responsible for the smallness of the overshoot. 

Yet we do not claim that our model is the final word 
concerning overshoot in the Sun, as several simplifiy- 
ing assumptions were made both in the model and in 
its numerical implementation. In fact, a definitely bet- 
ter closure for the non-local fluxes is already available 
(ICanutoU200llh and it has been app lied to A stars success- 
fully ([Kupka fc Montgomery! 12002). Further work is also 
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Fig. 3. Distribution of turbulent kinetic energy fc, convec- 
tive flux equivalent speed J, turbulent length scale I, and 
superadiabaticity AV in an (unrealistic) solution com- 
puted with a simple downgradient closure, for comparison. 



nee ded to inc l ude a consistent treatment of the anisotropy 
Icf. IPetrovavl Il992t . to clarify the dependence of the re- 
sults on parameter values and closure assumptions, and 
to extend the domain of calculation to the unstable layer. 
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